Info from Kim:

  • index from 1 to 6 stands for groups where: 1- HBC, 2-iOSN, 3-HBC*, 4-GBC, 5-mOSN, 6-all groups
  • names with “dataTFtranscel” include only TF genes
  • names with “data1500transcel” stand for data that have 1500 genes (mixing genes)
  • all results received by performing zinb1 algorithm with alpha(significance levels tests) \(\alpha=2\times\text{pnorm}(n^{0.2},\text{lower.tail=F})\) ; except filenames that have “alpha059515result”; here \(\alpha=2\times\text{pnorm}(n^{0.15},\text{lower.tail=F})\)
library(igraph)
## 
## Attaching package: 'igraph'
## The following objects are masked from 'package:stats':
## 
##     decompose, spectrum
## The following object is masked from 'package:base':
## 
##     union
library(rgl)
library(rafalib)
library(slingshot)
## Loading required package: princurve
library(msigdbr)
library(fgsea)
## Loading required package: Rcpp
library(knitr)
library(leiden)
sds <- readRDS("../data/finalTrajectory/sling.rds")
clDatta <- readRDS("../data/finalTrajectory/dattaCl_noResp_noMV.rds")
cols <- c("#E41A1C", "#377EB8", "#4DAF4A", "#984EA3", "#FF7F00", "#FFFF33", "#A65628", "#F781BF", "#999999")
view3d( theta = 10, phi = 10)
# Datta labels
plot3d(reducedDim(sds), aspect = 'iso', col=cols[clDatta], alpha=.6, xlab="UMAP1", ylab="UMAP2", zlab="UMAP3", box=FALSE, top=TRUE)
load("../data/list_TF_Kim.RData")

Transcription factor networks

Neuronal lineage

HBC* - GBC - iOSN - mOSN

## activated HBC
load("../data/associationNetworkKim/result_realdata/3-dataTFtranscelsalpha059515result.RData")
adjHBCAct <- adj.zinb1
dimnames(adjHBCAct) <- list(colnames(datatest), colnames(datatest))
adjHBCAct <- adjHBCAct[!rowSums(adjHBCAct) == 0,]
adjHBCAct <- adjHBCAct[,!colSums(adjHBCAct) == 0]
rm(adj.zinb1, datatest)

## GBC
load("../data/associationNetworkKim/result_realdata/4-dataTFtranscelsalpha059515result.RData")
adjGBC <- adj.zinb1
dimnames(adjGBC) <- list(colnames(datatest), colnames(datatest))
adjGBC <- adjGBC[!rowSums(adjGBC) == 0,]
adjGBC <- adjGBC[,!colSums(adjGBC) == 0]
rm(adj.zinb1, datatest)

## iOSN
load("../data/associationNetworkKim/result_realdata/2-dataTFtranscelsalpha059515result.RData")
adjiOSN <- adj.zinb1
dimnames(adjiOSN) <- list(colnames(datatest), colnames(datatest))
adjiOSN <- adjiOSN[!rowSums(adjiOSN) == 0,]
adjiOSN <- adjiOSN[,!colSums(adjiOSN) == 0]
rm(adj.zinb1, datatest)

## mOSN
load("../data/associationNetworkKim/result_realdata/5-dataTFtranscelsalpha059515result.RData")
adjmOSN <- adj.zinb1
dimnames(adjmOSN) <- list(colnames(datatest), colnames(datatest))
adjmOSN <- adjmOSN[!rowSums(adjmOSN) == 0,]
adjmOSN <- adjmOSN[,!colSums(adjmOSN) == 0]
rm(adj.zinb1, datatest)
# gHBC <- graph_from_adjacency_matrix(adjHBC,mode="undirected")
# cols <- c(brewer.pal(8, "Dark2"),brewer.pal(8, "Set2"))
# V(g)$color <- cols[cl$cluster[order(cl$X, decreasing=FALSE)]+2]
# lo <- layout_with_fr(gHBC)
# plot(g,vertex.size=3,vertex.label=NA,vertex.frame.color=V(g)$color,layout=lo,edge.width=0.6)
# legend("bottomleft",paste("comm",seq(1,6)),col=cols[3:8],pch=16, bty='n')

degreeDf <- data.frame(degree=1:52,
                       hbcAct=NA,
                       gbc=NA,
                       iOSN=NA,
                       mOSN=NA)
degreeHBCAct <- table(rowSums(adjHBCAct))
degreeDf$hbcAct[match(names(degreeHBCAct), degreeDf$degree)] <- degreeHBCAct
degreeGBC <- table(rowSums(adjGBC))
degreeDf$gbc[match(names(degreeGBC), degreeDf$degree)] <- degreeGBC
degreeiOSN <- table(rowSums(adjiOSN))
degreeDf$iOSN[match(names(degreeiOSN), degreeDf$degree)] <- degreeiOSN
degreemOSN <- table(rowSums(adjmOSN))
degreeDf$mOSN[match(names(degreemOSN), degreeDf$degree)] <- degreemOSN


mypar(mfrow=c(2,2))
barplot(degreeDf$hbcAct, names=degreeDf$degree, main="HBC*: Node-level degree", ylim=c(0,400))
barplot(degreeDf$gbc, names=degreeDf$degree, main="GBC: Node-level degree", ylim=c(0,400))
barplot(degreeDf$iOSN, names=degreeDf$degree, main="iOSN: Node-level degree", ylim=c(0,400))
barplot(degreeDf$mOSN, names=degreeDf$degree, main="mOSN: Node-level degree", ylim=c(0,400))

Extract 2-core of the network

extractCore <- function(A, degree = 3){
  converg <- FALSE
  old.nrow <- nrow(A)
  while(!converg){
    d <- colSums(A)
    to.keep <- which(d>=degree)
    if(old.nrow==length(to.keep)){
      converg <- TRUE
    }
    old.nrow <- length(to.keep)
    A <- A[to.keep,to.keep]
  }
  return(A)
}

adjHBCAct2 <- extractCore(adjHBCAct, degree = 2)
adjGBC2 <- extractCore(adjGBC, degree = 2)
adjiOSN2 <- extractCore(adjiOSN, degree = 2)
adjmOSN2 <- extractCore(adjmOSN, degree = 2)

mypar(mfrow=c(2,2))
gHBCAct2 <- graph_from_adjacency_matrix(adjHBCAct2,mode="undirected")
plot(gHBCAct2,vertex.size=3,vertex.label=NA,edge.width=0.6, main="HBC*")
gGBC2 <- graph_from_adjacency_matrix(adjGBC2,mode="undirected")
plot(gGBC2,vertex.size=3,vertex.label=NA,edge.width=0.6, main="GBC")
giOSN2 <- graph_from_adjacency_matrix(adjiOSN2,mode="undirected")
plot(gGBC2,vertex.size=3,vertex.label=NA,edge.width=0.6, main="iOSN")
gmOSN2 <- graph_from_adjacency_matrix(adjmOSN2,mode="undirected")
plot(gGBC2,vertex.size=3,vertex.label=NA,edge.width=0.6, main="mOSN")

saveRDS(list("hbcAct" = adjHBCAct2, 
             "GBC" = adjGBC2, 
             "iOSN" = adjiOSN2, 
             "mOSN" = adjmOSN2), file="../data/adjMatrices_2core.rds")

Community detection

library(RColorBrewer)
library(HCD)

gradualPlot <- function(cl, A){
  ncl <- length(unique(cl))
  for(cc in 1:ncl){
    id <- which(cl %in% (1:cc))
    curG <- graph_from_adjacency_matrix(A[id,id],mode="undirected")
      V(curG)$color <- cols[cl[names(V(curG))]]
      lo <- layout_with_fr(curG)
     plot(curG,
         vertex.size=3,
       vertex.label=NA,
       edge.width=0.6,
       vertex.frame.color=cols[cl[id]],
       color=V(curG)$color,
       layout = lo
       )
  }
}

## C5 category is according to gene ontology grouping: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC4707969/pdf/nihms-743907.pdf
geneSets <- msigdbr(species = "Mus musculus", category = "C5", subcategory = "BP")

gsea <- function(genes, background, geneSets, n=10, minSize=5, name=NULL){
  ### filter background to only include genes that we assessed.
  geneSets <- geneSets[geneSets$gene_symbol %in% background,]
  m_list <- geneSets %>% split(x = .$gene_symbol, f = .$gs_name)
  # gene set must have at least minSize genes in background.
  m_list <- m_list[unlist(lapply(m_list, length)) >= minSize]
  
  overlapPval <- unlist(lapply(m_list, function(gs){
    # genes in community and gene set
    inBoth <- sum(genes %in% gs)
    # genes in community and not in gene set
    inComOnly <- length(genes) - inBoth
    # genes in background and gene set
    inGsBack <- sum(background %in% gs)
    # genes in background and not in gene set
    outGsBack <- length(background) - inGsBack
    m <- matrix(c(inBoth, inComOnly,
           inGsBack, outGsBack),
           nrow =2, ncol=2, byrow=TRUE,
           dimnames = list(c("in community", "out community"),
                           c("in gene set", "out gene set")))
    fis <- fisher.test(m, alternative = "greater")
    pval <- fis$p.value
    return(pval)
  }))
  padj <- p.adjust(overlapPval, "fdr")
  oo <- order(overlapPval, decreasing=FALSE)
  res <- data.frame(geneSet = names(m_list)[oo[1:n]],
                    pval = overlapPval[oo[1:n]],
                    padj = padj[oo[1:n]],
                    row.names = NULL)
  kable(res, caption=name, label=name)
}

HBC*

#set.seed(2)
clHBCAct <- leiden(adjHBCAct2, resolution_parameter=.4, seed = 7)
table(clHBCAct)
## clHBCAct
##   1   2   3   4 
## 378 286 234  28
gHBCAct <- graph_from_adjacency_matrix(adjHBCAct2,mode="undirected")
cols <- brewer.pal(8, "Dark2")
V(gHBCAct)$color <- cols[clHBCAct[names(V(gHBCAct))]]
plot(gHBCAct,
     vertex.size=3,
     vertex.label=NA,
     #vertex.frame.color=V(gHBC)$color,
     color=V(gHBCAct)$color,
     edge.width=0.6)

comHBCAct <- make_clusters(gHBCAct, membership = clHBCAct)
plot(comHBCAct, gHBCAct,vertex.size=3, vertex.label=NA, edge.width=0.6)

## gradual plot
gradualPlot(clHBCAct, adjHBCAct2)

## gene set enrichment on HCD results
ncl <- length(unique(comHBCAct$membership))
kab <- list()
# kabNames <- c("Krt / cell cycle / stimulus",
#                 "Chromatin",
#                 "Epigenetics, cell cycle",
#                 "Cell cycle")
kabNames <- c("Cell cycle, DNA damage",
                "Epigenetics",
                "(Epithelial) cell differentiation",
                "Cell cycle, DNA replication",
              "DNA replication")
for(kk in 1:ncl){
  genes <- rownames(adjHBCAct2)[which(comHBCAct$membership == kk)]
  kab[[kk]] <- gsea(genes = genes,
       background = TF.list,
       geneSets = geneSets,
       name = kabNames[kk])
}
kab
## [[1]]
## 
## 
## Table: Cell cycle, DNA damage
## 
## geneSet                                              pval        padj
## ---------------------------------------------  ----------  ----------
## GO_CELLULAR_RESPONSE_TO_DNA_DAMAGE_STIMULUS     0.0010309   0.6004173
## GO_REGULATION_OF_CELL_CYCLE                     0.0023338   0.6004173
## GO_CARBOHYDRATE_HOMEOSTASIS                     0.0023584   0.6004173
## GO_RESPONSE_TO_CARBOHYDRATE                     0.0029978   0.6004173
## GO_PROCESS_UTILIZING_AUTOPHAGIC_MECHANISM       0.0033514   0.6004173
## GO_REGULATION_OF_MITOTIC_CELL_CYCLE             0.0034382   0.6004173
## GO_NEGATIVE_REGULATION_OF_CELL_CYCLE_PROCESS    0.0037266   0.6004173
## GO_REGULATION_OF_AUTOPHAGY                      0.0041650   0.6004173
## GO_DNA_REPAIR                                   0.0044003   0.6004173
## GO_REGULATION_OF_CELL_CYCLE_PHASE_TRANSITION    0.0045432   0.6004173
## 
## [[2]]
## 
## 
## Table: Epigenetics
## 
## geneSet                                                         pval   padj
## --------------------------------------------------------  ----------  -----
## GO_PROTEIN_DEMETHYLATION                                   0.0137397      1
## GO_DEMETHYLATION                                           0.0196971      1
## GO_REGULATION_OF_DNA_REPLICATION                           0.0275424      1
## GO_FOAM_CELL_DIFFERENTIATION                               0.0278109      1
## GO_HISTONE_H3_K4_DEMETHYLATION                             0.0385626      1
## GO_HISTONE_H3_K4_DEMETHYLATION_TRIMETHYL_H3_K4_SPECIFIC    0.0385626      1
## GO_HISTONE_H3_K4_METHYLATION                               0.0479833      1
## GO_RIBONUCLEOPROTEIN_COMPLEX_SUBUNIT_ORGANIZATION          0.0550234      1
## GO_REGULATION_OF_CHOLESTEROL_BIOSYNTHETIC_PROCESS          0.0708770      1
## GO_STEROL_BIOSYNTHETIC_PROCESS                             0.0708770      1
## 
## [[3]]
## 
## 
## Table: (Epithelial) cell differentiation
## 
## geneSet                                                 pval        padj
## ------------------------------------------------  ----------  ----------
## GO_SCHWANN_CELL_DIFFERENTIATION                    0.0037526   0.9447078
## GO_REGULATION_OF_CELL_POPULATION_PROLIFERATION     0.0051764   0.9447078
## GO_EPIDERMIS_DEVELOPMENT                           0.0064884   0.9447078
## GO_POSITIVE_REGULATION_OF_CELL_DEATH               0.0072210   0.9447078
## GO_PERIPHERAL_NERVOUS_SYSTEM_DEVELOPMENT           0.0106454   0.9447078
## GO_RAS_PROTEIN_SIGNAL_TRANSDUCTION                 0.0117230   0.9447078
## GO_EPITHELIUM_DEVELOPMENT                          0.0134552   0.9447078
## GO_NEGATIVE_REGULATION_OF_BMP_SIGNALING_PATHWAY    0.0137142   0.9447078
## GO_REGULATION_OF_CELL_CYCLE                        0.0159675   0.9447078
## GO_APOPTOTIC_PROCESS                               0.0161161   0.9447078
## 
## [[4]]
## 
## 
## Table: Cell cycle, DNA replication
## 
## geneSet                                              pval        padj
## ---------------------------------------------  ----------  ----------
## GO_DNA_REPLICATION_INITIATION                   0.0000000   0.0001121
## GO_DNA_DEPENDENT_DNA_REPLICATION                0.0000081   0.0075936
## GO_CELL_CYCLE_G1_S_PHASE_TRANSITION             0.0000101   0.0075936
## GO_DNA_UNWINDING_INVOLVED_IN_DNA_REPLICATION    0.0000158   0.0089155
## GO_DNA_GEOMETRIC_CHANGE                         0.0001095   0.0494663
## GO_CELL_CYCLE_PHASE_TRANSITION                  0.0002037   0.0766622
## GO_DNA_REPLICATION                              0.0003546   0.1041845
## GO_MITOTIC_CELL_CYCLE                           0.0003691   0.1041845
## GO_DNA_METABOLIC_PROCESS                        0.0009372   0.2351229
## GO_DNA_CONFORMATION_CHANGE                      0.0013729   0.3100113

GBC

set.seed(3)
clGBC <- leiden(adjGBC2, resolution_parameter=.4, seed=7)
table(clGBC)
## clGBC
##   1   2   3   4 
## 331 279 274 107
gGBC <- graph_from_adjacency_matrix(adjGBC2,mode="undirected")
cols <- brewer.pal(8, "Dark2")
V(gGBC)$color <- cols[clGBC[names(V(gGBC))]]
plot(gGBC,
     vertex.size=3,
     vertex.label=NA,
     #vertex.frame.color=V(gHBC)$color,
     color=V(gGBC)$color,
     edge.width=0.6)

comGBC <- make_clusters(gGBC, membership = clGBC)
plot(comGBC, gGBC,vertex.size=3, vertex.label=NA, edge.width=0.6)

## gradual plot
gradualPlot(clGBC, adjGBC2)

## gene set enrichment 
ncl <- length(unique(comGBC$membership))
kab <- list()
# kabNames <- c("Chromatin, neuron/cell development",
#                 "P53, DNA replication, cell cycle",
#                 "Cell differentiation / signaling",
#                 "Notch, signaling",
#               "DNA repair")
kabNames <- c("P53, DNA replication",
                "Notch signaling",
                "Cell differentiation / signaling",
                "RNA splicing, DNA repair, epigenetics",
                "NIK/NF Kappa-B signaling")
for(kk in 1:ncl){
  genes <- rownames(adjGBC2)[which(comGBC$membership == kk)]
  print(gsea(genes = genes,
       background = TF.list,
       geneSets = geneSets,
       name = kabNames[kk]))
}
## 
## 
## Table: P53, DNA replication
## 
## geneSet                                                             pval   padj
## ------------------------------------------------------------  ----------  -----
## GO_REGULATION_OF_SIGNAL_TRANSDUCTION_BY_P53_CLASS_MEDIATOR     0.0073943      1
## GO_DNA_REPLICATION_INITIATION                                  0.0088812      1
## GO_HISTONE_PHOSPHORYLATION                                     0.0088812      1
## GO_SIGNAL_TRANSDUCTION_BY_P53_CLASS_MEDIATOR                   0.0111609      1
## GO_DNA_REPLICATION                                             0.0173178      1
## GO_SNRNA_TRANSCRIPTION                                         0.0176173      1
## GO_NCRNA_TRANSCRIPTION                                         0.0184906      1
## GO_TRANSCRIPTION_ELONGATION_FROM_RNA_POLYMERASE_II_PROMOTER    0.0212404      1
## GO_POSITIVE_REGULATION_OF_CHROMOSOME_ORGANIZATION              0.0277202      1
## GO_DNA_TEMPLATED_TRANSCRIPTION_ELONGATION                      0.0301290      1
## 
## 
## Table: Notch signaling
## 
## geneSet                                                pval        padj
## -----------------------------------------------  ----------  ----------
## GO_REGULATION_OF_CELL_POPULATION_PROLIFERATION    0.0003936   0.3441872
## GO_REGULATION_OF_NOTCH_SIGNALING_PATHWAY          0.0004437   0.3441872
## GO_RESPONSE_TO_NITROGEN_COMPOUND                  0.0005204   0.3441872
## GO_NEGATIVE_REGULATION_OF_SIGNALING               0.0007461   0.3441872
## GO_CELLULAR_RESPONSE_TO_EXTERNAL_STIMULUS         0.0010464   0.3441872
## GO_RESPONSE_TO_PEPTIDE                            0.0013035   0.3441872
## GO_RESPONSE_TO_STARVATION                         0.0016670   0.3441872
## GO_CELLULAR_RESPONSE_TO_STARVATION                0.0016766   0.3441872
## GO_CARDIOVASCULAR_SYSTEM_DEVELOPMENT              0.0017096   0.3441872
## GO_NOTCH_SIGNALING_PATHWAY                        0.0020545   0.3441872
## 
## 
## Table: Cell differentiation / signaling
## 
## geneSet                                                   pval   padj
## --------------------------------------------------  ----------  -----
## GO_RNA_SPLICING_VIA_TRANSESTERIFICATION_REACTIONS    0.0055550      1
## GO_REGULATION_OF_GENE_SILENCING                      0.0070297      1
## GO_RNA_SPLICING                                      0.0077460      1
## GO_CELLULAR_PROTEIN_CONTAINING_COMPLEX_ASSEMBLY      0.0084200      1
## GO_DNA_REPAIR                                        0.0088367      1
## GO_MRNA_METABOLIC_PROCESS                            0.0104386      1
## GO_TRANSCRIPTION_PREINITIATION_COMPLEX_ASSEMBLY      0.0106774      1
## GO_REGULATION_OF_GENE_EXPRESSION_EPIGENETIC          0.0139632      1
## GO_REGULATION_OF_STEM_CELL_DIFFERENTIATION           0.0142288      1
## GO_REGULATION_OF_DNA_BINDING                         0.0146403      1
## 
## 
## Table: RNA splicing, DNA repair, epigenetics
## 
## geneSet                                                              pval   padj
## -------------------------------------------------------------  ----------  -----
## GO_REGULATION_OF_NIK_NF_KAPPAB_SIGNALING                        0.0034495      1
## GO_ORGANELLE_LOCALIZATION                                       0.0051168      1
## GO_NIK_NF_KAPPAB_SIGNALING                                      0.0164556      1
## GO_ESTABLISHMENT_OF_ORGANELLE_LOCALIZATION                      0.0206255      1
## GO_MICROTUBULE_CYTOSKELETON_ORGANIZATION_INVOLVED_IN_MITOSIS    0.0227626      1
## GO_POSITIVE_REGULATION_OF_NIK_NF_KAPPAB_SIGNALING               0.0298471      1
## GO_SPINDLE_ASSEMBLY                                             0.0298471      1
## GO_HISTONE_H3_ACETYLATION                                       0.0309471      1
## GO_HISTONE_H4_ACETYLATION                                       0.0355192      1
## GO_CELLULAR_RESPONSE_TO_CALCIUM_ION                             0.0379565      1
kab
## list()

iOSN

set.seed(4)
cliOSN <- leiden(adjiOSN2, resolution_parameter=.3, seed=7)
table(cliOSN)
## cliOSN
##   1   2   3   4   5   6   7 
## 235  83  77  27  21  16   6
giOSN <- graph_from_adjacency_matrix(adjiOSN2,mode="undirected")
cols <- brewer.pal(8, "Dark2")
V(giOSN)$color <- cols[cliOSN[names(V(giOSN))]]
plot(giOSN,
     vertex.size=3,
     vertex.label=NA,
     color=V(giOSN)$color,
     edge.width=0.6)

# igraph plot
comiOSN <- make_clusters(giOSN, membership = cliOSN)
plot(comiOSN, giOSN,vertex.size=3, vertex.label=NA, edge.width=0.6)

## gradual plot
gradualPlot(cliOSN, adjiOSN2)

## gene set enrichment on HCD results
ncl <- length(unique(comiOSN$membership))
kab <- list()
kabNames <- c("Neurogenesis",
              "Defense (virus) response",
              "P53, membrane",
              "Axon extension, cytokinesis, epigenetic gene expression regulation",
              "Cytokinesis, DNA replication",
              "Interleukins, cell cycle, protein localization",
              "Stem cell population maintenance, transcription")
for(kk in 1:ncl){
  genes <- rownames(adjiOSN2)[which(comiOSN$membership == kk)]
  print(gsea(genes = genes,
       background = TF.list,
       geneSets = geneSets,
       name = kabNames[kk]))
}
## 
## 
## Table: Neurogenesis
## 
## geneSet                                                         pval        padj
## --------------------------------------------------------  ----------  ----------
## GO_TELENCEPHALON_DEVELOPMENT                               0.0002684   0.3249968
## GO_NEUROGENESIS                                            0.0002879   0.3249968
## GO_NEURON_DIFFERENTIATION                                  0.0008345   0.3692460
## GO_REGULATION_OF_NERVOUS_SYSTEM_DEVELOPMENT                0.0008941   0.3692460
## GO_EMBRYO_DEVELOPMENT                                      0.0008958   0.3692460
## GO_CENTRAL_NERVOUS_SYSTEM_DEVELOPMENT                      0.0011003   0.3692460
## GO_POSITIVE_REGULATION_OF_CELL_POPULATION_PROLIFERATION    0.0011447   0.3692460
## GO_FOREBRAIN_DEVELOPMENT                                   0.0014588   0.3837962
## GO_NEGATIVE_REGULATION_OF_NEURON_DIFFERENTIATION           0.0015983   0.3837962
## GO_REGIONALIZATION                                         0.0019554   0.3837962
## 
## 
## Table: Defense (virus) response
## 
## geneSet                                                                    pval   padj
## -------------------------------------------------------------------  ----------  -----
## GO_DEFENSE_RESPONSE_TO_VIRUS                                          0.0032635      1
## GO_RESPONSE_TO_VIRUS                                                  0.0038295      1
## GO_RESPONSE_TO_BIOTIC_STIMULUS                                        0.0056459      1
## GO_NEGATIVE_REGULATION_OF_JNK_CASCADE                                 0.0118638      1
## GO_REGULATION_OF_MRNA_CATABOLIC_PROCESS                               0.0119217      1
## GO_REGULATION_OF_STRESS_ACTIVATED_PROTEIN_KINASE_SIGNALING_CASCADE    0.0119217      1
## GO_RESPONSE_TO_TYPE_I_INTERFERON                                      0.0135956      1
## GO_CELLULAR_HOMEOSTASIS                                               0.0159622      1
## GO_DEFENSE_RESPONSE_TO_OTHER_ORGANISM                                 0.0176323      1
## GO_RESPONSE_TO_CYTOKINE                                               0.0185354      1
## 
## 
## Table: P53, membrane
## 
## geneSet                                                                                          pval   padj
## -----------------------------------------------------------------------------------------  ----------  -----
## GO_MITOCHONDRIAL_MEMBRANE_ORGANIZATION                                                      0.0209248      1
## GO_PROTEIN_MONOUBIQUITINATION                                                               0.0209248      1
## GO_REGULATION_OF_MEMBRANE_PERMEABILITY                                                      0.0209248      1
## GO_INTRINSIC_APOPTOTIC_SIGNALING_PATHWAY_BY_P53_CLASS_MEDIATOR                              0.0242866      1
## GO_INTRINSIC_APOPTOTIC_SIGNALING_PATHWAY_IN_RESPONSE_TO_DNA_DAMAGE_BY_P53_CLASS_MEDIATOR    0.0310897      1
## GO_MAINTENANCE_OF_PROTEIN_LOCATION                                                          0.0310897      1
## GO_REGULATION_OF_GENE_SILENCING                                                             0.0312073      1
## GO_MRNA_TRANSCRIPTION                                                                       0.0369826      1
## GO_REGULATION_OF_CHROMATIN_SILENCING                                                        0.0369826      1
## GO_BASE_EXCISION_REPAIR                                                                     0.0400952      1
## 
## 
## Table: Axon extension, cytokinesis, epigenetic gene expression regulation
## 
## geneSet                                                  pval   padj
## -------------------------------------------------  ----------  -----
## GO_REGULATION_OF_VESICLE_MEDIATED_TRANSPORT         0.0718594      1
## GO_DNA_TEMPLATED_TRANSCRIPTION_TERMINATION          0.0870021      1
## GO_REGULATION_OF_GENE_EXPRESSION_EPIGENETIC         0.0879146      1
## GO_ATF6_MEDIATED_UNFOLDED_PROTEIN_RESPONSE          0.0989998      1
## GO_AXON_EXTENSION                                   0.0989998      1
## GO_CHOLESTEROL_EFFLUX                               0.0989998      1
## GO_CHOLESTEROL_STORAGE                              0.0989998      1
## GO_CYTOKINESIS                                      0.0989998      1
## GO_MITOTIC_SPINDLE_ASSEMBLY                         0.0989998      1
## GO_NEGATIVE_REGULATION_OF_INNATE_IMMUNE_RESPONSE    0.0989998      1
## 
## 
## Table: Cytokinesis, DNA replication
## 
## geneSet                                                                    pval   padj
## -------------------------------------------------------------------  ----------  -----
## GO_CYTOKINESIS                                                        0.0034643      1
## GO_POSITIVE_REGULATION_OF_DNA_DEPENDENT_DNA_REPLICATION               0.0034643      1
## GO_REGULATION_OF_CYTOKINESIS                                          0.0034643      1
## GO_TROPHOBLAST_GIANT_CELL_DIFFERENTIATION                             0.0034643      1
## GO_EMBRYO_DEVELOPMENT_ENDING_IN_BIRTH_OR_EGG_HATCHING                 0.0062909      1
## GO_CELL_DIFFERENTIATION_INVOLVED_IN_EMBRYONIC_PLACENTA_DEVELOPMENT    0.0072449      1
## GO_POSITIVE_REGULATION_OF_DNA_REPLICATION                             0.0072449      1
## GO_POSITIVE_REGULATION_OF_CELL_CYCLE_ARREST                           0.0087899      1
## GO_EMBRYONIC_PLACENTA_DEVELOPMENT                                     0.0095574      1
## GO_CELL_CYCLE_DNA_REPLICATION                                         0.0104550      1
## 
## 
## Table: Interleukins, cell cycle, protein localization
## 
## geneSet                                                                 pval   padj
## ----------------------------------------------------------------  ----------  -----
## GO_RESPONSE_TO_INTERLEUKIN_4                                       0.0072146      1
## GO_G0_TO_G1_TRANSITION                                             0.0168525      1
## GO_CIRCADIAN_RHYTHM                                                0.0284399      1
## GO_POSITIVE_REGULATION_OF_CELLULAR_PROTEIN_LOCALIZATION            0.0341263      1
## GO_DNA_INTEGRITY_CHECKPOINT                                        0.0482564      1
## GO_RHYTHMIC_PROCESS                                                0.0508153      1
## GO_POSITIVE_REGULATION_OF_ESTABLISHMENT_OF_PROTEIN_LOCALIZATION    0.0560070      1
## GO_ESTABLISHMENT_OF_PROTEIN_LOCALIZATION_TO_ORGANELLE              0.0586845      1
## GO_CHROMATIN_SILENCING_AT_TELOMERE                                 0.0601135      1
## GO_HISTONE_H3_K9_ACETYLATION                                       0.0601135      1
## 
## 
## Table: Stem cell population maintenance, transcription
## 
## geneSet                                                                     pval   padj
## --------------------------------------------------------------------  ----------  -----
## GO_REGULATION_OF_STEM_CELL_POPULATION_MAINTENANCE                      0.0306392      1
## GO_REGULATION_OF_DNA_TEMPLATED_TRANSCRIPTION_ELONGATION                0.0456629      1
## GO_REGULATION_OF_PROTEIN_COMPLEX_ASSEMBLY                              0.1352331      1
## GO_REGULATION_OF_CELLULAR_AMIDE_METABOLIC_PROCESS                      0.1386647      1
## GO_PEPTIDE_BIOSYNTHETIC_PROCESS                                        0.1522779      1
## GO_AMIDE_BIOSYNTHETIC_PROCESS                                          0.1623697      1
## GO_DNA_TEMPLATED_TRANSCRIPTION_ELONGATION                              0.1623697      1
## GO_POSITIVE_REGULATION_OF_DNA_BINDING_TRANSCRIPTION_FACTOR_ACTIVITY    0.1623697      1
## GO_CELLULAR_AMIDE_METABOLIC_PROCESS                                    0.1756696      1
## GO_MAINTENANCE_OF_CELL_NUMBER                                          0.1952891      1
kab
## list()

mOSN

set.seed(2)
clmOSN <- leiden(adjmOSN2, resolution_parameter=.4, seed=7)
table(clmOSN)
## clmOSN
##   1   2   3   4   5 
## 329 182 139  81  41
gmOSN <- graph_from_adjacency_matrix(adjmOSN2,mode="undirected")
cols <- brewer.pal(8, "Dark2")
V(gmOSN)$color <- cols[clmOSN[names(V(gmOSN))]]
plot(gmOSN,
     vertex.size=3,
     vertex.label=NA,
     color=V(gmOSN)$color,
     edge.width=0.6)

# igraph plot
commOSN <- make_clusters(gmOSN, membership = clmOSN)
plot(commOSN, gmOSN,vertex.size=3, vertex.label=NA, edge.width=0.6)

## gradual plot
gradualPlot(clmOSN, adjmOSN2)

## gene set enrichment on HCD results
ncl <- length(unique(commOSN$membership))
kab <- list()
kabNames <- c("Chromatin",
              "Chromatin, transcription",
              "TGBF-Beta pathway, neuron death",
              "Viral processes",
              "Cell projection, others")
for(kk in 1:ncl){
  genes <- rownames(adjmOSN2)[which(commOSN$membership == kk)]
  print(gsea(genes = genes,
       background = TF.list,
       geneSets = geneSets,
       name = kabNames[kk]))
}
## 
## 
## Table: Chromatin
## 
## geneSet                                                   pval   padj
## --------------------------------------------------  ----------  -----
## GO_COVALENT_CHROMATIN_MODIFICATION                   0.0076200      1
## GO_CHROMATIN_ORGANIZATION                            0.0088562      1
## GO_PROTEIN_ACYLATION                                 0.0090561      1
## GO_CHROMOSOME_ORGANIZATION                           0.0124554      1
## GO_RNA_SPLICING_VIA_TRANSESTERIFICATION_REACTIONS    0.0151444      1
## GO_RRNA_TRANSCRIPTION                                0.0160161      1
## GO_POSITIVE_REGULATION_OF_DEFENSE_RESPONSE           0.0176331      1
## GO_PROTEIN_ACETYLATION                               0.0179776      1
## GO_PEPTIDYL_LYSINE_ACETYLATION                       0.0196668      1
## GO_REGULATION_OF_RNA_SPLICING                        0.0210600      1
## 
## 
## Table: Chromatin, transcription
## 
## geneSet                                                                        pval   padj
## -----------------------------------------------------------------------  ----------  -----
## GO_ATP_DEPENDENT_CHROMATIN_REMODELING                                     0.0099950      1
## GO_DNA_TEMPLATED_TRANSCRIPTION_ELONGATION                                 0.0113455      1
## GO_PROTEIN_DNA_COMPLEX_SUBUNIT_ORGANIZATION                               0.0202977      1
## GO_FERTILIZATION                                                          0.0218270      1
## GO_RESPONSE_TO_EXTRACELLULAR_STIMULUS                                     0.0219166      1
## GO_RESPONSE_TO_OSMOTIC_STRESS                                             0.0300984      1
## GO_CHROMATIN_ORGANIZATION                                                 0.0313435      1
## GO_INSULIN_RECEPTOR_SIGNALING_PATHWAY                                     0.0338212      1
## GO_ORGANOPHOSPHATE_BIOSYNTHETIC_PROCESS                                   0.0338212      1
## GO_CELLULAR_PROCESS_INVOLVED_IN_REPRODUCTION_IN_MULTICELLULAR_ORGANISM    0.0358365      1
## 
## 
## Table: TGBF-Beta pathway, neuron death
## 
## geneSet                                                                             pval        padj
## ----------------------------------------------------------------------------  ----------  ----------
## GO_TRANSFORMING_GROWTH_FACTOR_BETA_RECEPTOR_SIGNALING_PATHWAY                  0.0001683   0.3800204
## GO_CELLULAR_RESPONSE_TO_INORGANIC_SUBSTANCE                                    0.0005379   0.6072730
## GO_RESPONSE_TO_TRANSFORMING_GROWTH_FACTOR_BETA                                 0.0020270   0.8644903
## GO_REGULATION_OF_DNA_BINDING                                                   0.0021189   0.8644903
## GO_NEURON_DEATH                                                                0.0027505   0.8644903
## GO_REGULATION_OF_DNA_TEMPLATED_TRANSCRIPTION_IN_RESPONSE_TO_STRESS             0.0027646   0.8644903
## GO_TRANSMEMBRANE_RECEPTOR_PROTEIN_SERINE_THREONINE_KINASE_SIGNALING_PATHWAY    0.0028099   0.8644903
## GO_POSITIVE_REGULATION_OF_NEURON_DEATH                                         0.0030629   0.8644903
## GO_CELLULAR_RESPONSE_TO_CALCIUM_ION                                            0.0038982   0.9780268
## GO_RESPONSE_TO_METAL_ION                                                       0.0049149   1.0000000
## 
## 
## Table: Viral processes
## 
## geneSet                                                                                           pval   padj
## ------------------------------------------------------------------------------------------  ----------  -----
## GO_VIRAL_LIFE_CYCLE                                                                          0.0010567      1
## GO_INTERSPECIES_INTERACTION_BETWEEN_ORGANISMS                                                0.0016139      1
## GO_PROTEIN_MONOUBIQUITINATION                                                                0.0039224      1
## GO_REGULATION_OF_VIRAL_LIFE_CYCLE                                                            0.0051478      1
## GO_PROTEIN_MODIFICATION_BY_SMALL_PROTEIN_CONJUGATION                                         0.0060424      1
## GO_REGULATION_OF_NUCLEAR_TRANSCRIBED_MRNA_CATABOLIC_PROCESS_DEADENYLATION_DEPENDENT_DECAY    0.0080769      1
## GO_REGULATION_OF_SYMBIOSIS_ENCOMPASSING_MUTUALISM_THROUGH_PARASITISM                         0.0103493      1
## GO_PROTEIN_SUMOYLATION                                                                       0.0125485      1
## GO_POSITIVE_REGULATION_OF_VIRAL_PROCESS                                                      0.0166502      1
## GO_CELLULAR_MACROMOLECULE_CATABOLIC_PROCESS                                                  0.0180804      1
## 
## 
## Table: Cell projection, others
## 
## geneSet                                                            pval   padj
## -----------------------------------------------------------  ----------  -----
## GO_CELL_PROJECTION_ASSEMBLY                                   0.0135173      1
## GO_REGULATION_OF_CIRCADIAN_RHYTHM                             0.0137572      1
## GO_NEGATIVE_REGULATION_OF_SMOOTHENED_SIGNALING_PATHWAY        0.0165868      1
## GO_ENDOSOMAL_TRANSPORT                                        0.0209788      1
## GO_POSITIVE_REGULATION_OF_CELL_CYCLE_G2_M_PHASE_TRANSITION    0.0209788      1
## GO_SPINAL_CORD_PATTERNING                                     0.0209788      1
## GO_PROTEIN_DEPHOSPHORYLATION                                  0.0257976      1
## GO_PROTEIN_LOCALIZATION_TO_MEMBRANE                           0.0437886      1
## GO_CELL_CYCLE_G2_M_PHASE_TRANSITION                           0.0474675      1
## GO_MITOTIC_SISTER_CHROMATID_SEGREGATION                       0.0488742      1
kab
## list()